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SOME NEW NYSTROH INTEGRATORS 

By William M. Lear 
TRW 

1 . INTRODUCTION 

Nystrbm integrators are self-starting integrators used to integrate the 
second order, vector, differential equations x = f(t,j<) and x. = £(t,)(,^). 

Honuki, Reference 1, in a very extensive study has found the Nystrbm type in- 
tegrators to be greatly superior to most Runge-Kutta type integrators in terms 
of speed and accuracy when integrating ballistic missile trajectories. 

Nystrbm integrator parameters are defined by a set of nonlinear constraint 
equations. Frequently there are more parameters than there are equations. When 
this is the case, the investigator is free to choose some additional constraint 
equations. lAdditional equations are included for higher order (improved accuracy) 
integration when x^ = £(t). 

Two of the lower order integrators given here are not new, and are due to 
Nystrom, Reference 2. They are given so that the reader will have a single source 
for all orders of Nystrom type integrators. 



2. TWO FUNCTION INTEGRATION OF x = f(t,x) 

This integrator evaluates f twice and hence the terminology "two function 
integration". It provides third order integration of x and x. If x = f(t), 
the integration will still only be third order for £ and The integration 
parameters are due to Nystrora, Reference 2, and are the solution of six con- 
straint equations in six unknowns, which are given in Reference 1. 

J<1 = x„) 

- ATf(t^ + 62AT', Xj, + * a^aTJcp ■ 

2in+l " + O(aT^) 

4i+l ^ ^ ^1-N * O(aT^) 

Where AT = t^+l'^n integration coefficients are given by 

62 = 2/3 
a^ = 2/9 
= 1/4 
a 2 = 1/4 
= 1/4 
= 3/4 
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3. THREE FUNCTION INTEGRATION OF x = f (t.x) 

This integrator provides fourth order integration for both x and x. 
However, if x = f(t), the integrator will be fifth order for both and x. 
The integration parameters are due to Lear. 

= in + AT^ + AT(a^k^ + a2Ji2 + a^kj) 

><n+l = 2^n ^1-1 ^ hh ^3-3 

The integration coefficients are given by 


62= .6 -^JT6 = .35505 10257 
6^ = .6 +yf706 = .84494 89743 
a^ = .21-.6 V^= .06303 06154 
b^ = {.15+4 >4^)725 = .04519 18359 
b2 = (5.1+11 V^)/25 = .31177 75487 
= 1/9 = .11111 mil 

02 = (7+20 V^)/36 = .33052 72081 

03 = (7-20 VTo6)/36 = .05836 16809 
= 1/9 = .11111 11111 

p2 = (8+5 V^)/18 = .51248 58262 
33 = (8-5 V^)/18 = .37640 30627 
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4. FOUR FUNCTION INTEGRATION OF x = fCt,xl 

This integrator provides fifth order integration for x and However, 
if 5^ = f(t), the integrator will be seventh order for both x and x. The 
integration parameters are due to Lear. 

= aTf(tj^ + 52AT, 3 -|AT!i]) 

J<3 = ATf[t^ + ^^AT, + 53AT^ + AT(b^J^^ + b2J^)] 

^ AT{c-jJc-j + ^2]^ + ^3113)1 

2^+1 ^ 2in ■*■ + AT(o^k;j + 02 ]^ + a^jc^ + a^j^) 

4+1 " 4 ^1-1 ^ hh ^ ^3^3 ^Ah 

. The integration coefficients were obtained by numerically solving 17 non- 
linear equations in 17 unknowns. The solution is 

«2 = .21234 05385 

63 = .59053 31358 

64 = .91141 20406 
a^ = .022544 25214 
fa^ = -.00114 39805 
b2 = .17550 86728 
c^ = .11715 41673 
C2 = .13937 54710 
C3 = .15880 63156 



== .06250 00001 
«2 = .25901 73402 
Oj “ .15895 23623 
« .01953 02974 
B-, = .06250 00001 
Bg = .32884 43202 
B 3 = .38819 34687 
B 4 = .22046 22110 

Note that probably both and B-j are exactly .0625 = 1/16. Computer 
roundoff error probably prevented the correct solution. However, the above 
set of coefficients should be used as shown, since they satisfy many of the 
defining constraint equations exactly. 
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5. FIVE FUNCTION INTEGRATION OF x = l(t,x-) 

This integrator provides sixth order integration for x and x. There 
are 22 constraint equations in 24 unknowns, thus two more constraint equations 
were chosen so as to cause eighth order integration of f,(t). Unfortunately, 
the numerical solutions for the new constraint equations were unsatisfactory. 
Thus the integration constants shown here are only good to sixth order for 
integrating _x = f(t). 


Jil - 

kg = ATfCtj, + ^ a.j4tk^) 

+ A^A'T^ + ATi(.b-jjc-j +‘>b2JCg)3 


-2 
^3 

^ = ATT[.t^ + S^AT. 4 + 64AT4 
= ATf[t|^ + 5gAl, 

= >4, + AT4, + A'T-(a^k^ + agkg + 031(3 “4^ “5^^ 


J<3 = .*ATf[t^ + 53AT, 4 

r, + A^ATXy^ + AT;(C'jJCj + C gkg ^3^3 )] 

[, 4 + + dgkg + d^3 + d^4)0 


Vl 

4,+l 


The integration coefficients, due to Lear, are 


c, = 


cti - 


1 

,2 

^3 " 3 

s - 2 
®4 " 3 


'5 

= 1 

1 

8 

^1 

bg = 0 




_ 1 
“ 9 

Cg = 0 

c =1 
*^3 9 




= 0 

H - 8 

^2 " ■ 11 

^^3 “ 11 


'4 

_ 9 
' 22 

_ 11 

_ 4 

- 9 



_ 9 

120 

“2 ■ ■ W 

“3 " 20 


“4 

“ m 

_ 11 
■ 120 

®2 15 

^3 = ^4 = 

27 

40 

^5 

_ 11 
120 


~ 191 


“5 " ° 
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Monuki's integration coefficients (Reference 1) are 


^2 = 


.3 



«4 = 2/3 65 = 1 


= .045 = .18 b 2 = 0 

c^ = .13671 69639 
C 2 = .08047 55373 
c^ = .00502 97211 
= .00740 74074 
d 2 = .49023 56902 
d 3 = -.57037 03704 
= .57272 72727 
a., = .'08796 29630 
tt 2 = .33670 03367 
03 = -.23148 14815 
= .30681 81818 
= 0 

= .08796 29630 
^2 = .48100 04810 
$3 = -.57870 37037 
8 ^ = .92045 45455 
gg = .08928 57143 
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6 . SIX FUNCTION INTEGRATION OF x ^ f{t,x) 

This Integrator provides seventh order integration for _x and k. There 
were 33 constraint equations in 32 unknowns, impossible to solve at first 
glance. However, Monuki, Reference 1, using a brilliant set of assumptions, 
was able to show how to obtain a solution. 



The integration coefficients, due to Monuki”, are 

62 = .10654 17886 

63 = .21308 35772 
6 ^ = .59267 23008 
65 = .916 

= .972 

a., = .56755 76359*10'^ 
b^ = .07567 43515*10"^ 
b 2 = .15134 87029- 10'^ 
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= .14003 61674 
C2 = -.25447 80570 
C3 = .29007 21177 
= -1.02164 36141 
d2 = 2.65397 01073 
dg = -1.48615 90950 
d^ = .27336 06017 
= -20.40832 94915 

62 = 50.31431 81086 

63 = -32.30441 78724 
= 2.94949 60939 

63 = -.07867 48385 

= .06271 70177 . 
a2 - 0 

03 = .25968 74616 
= .15875’ 55586 
og = .01912 37845 
= -.00028 38224 
= .06271 70177 
$2 = 0 

^3 = .33000 64074 

64 = .38974 89881 
gg = .22766 41014 
pg = -.01013 65146 
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7. TWO FUNCTION INTEGRATION OF x = I(t,x,x) 

This integrator provides third order integration for x^ and second order 
for However, if x = £{t,x) or if }( = f.(t) then the integration will be 
third order for both x and x. The integration parameters are due to Nystrom, 
Reference 2, and are the solution of seven constraint equations in seven 
unknowns. 


i, “ x„. 4,) 

^ ^ ^ a^ATJc-j, ^ 

^+i " h'*' ^1-1 ®2-2 

The integration coefficients are given by 


62 = 2/3 
. = 2/9 

a-| = 1/4 
02 = 1/4 
a^ = 2/3 
= 1/4 
= 3/4 


0EIGWM.PAGEK 
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8. THREE FUNCTION INTEGRATION OF x = f (t,x,x) 

This integrator provides fourth order integration for x and third order 
for 3(. If = f.(t,)() the integrator will be fourth order for both and x. 

If j( = £(t) the integrator will be fifth order for both x and x. The inte- 
gration constants are due to Lear. 

= ATf{t^,2^,^) 

k2 = ATf(t|^ + SgAT, x„ + a^ATJip x^ + 

J<3 = ATf[t^ + 63A T, )^ + «3AT^ + AT(b^k^ + * ^iil ^2^^ 

>^+l = + AT^ + AT{a^k^ + ct2Jl2 + ^ 3 ^ 3 ) 

+ ^2-2 ^ ^3^3 

The integration coefficients are given by 

62 = .6- VToe = .35505 10257 

63 = .6+ypm= .84494 89743 

a-, = .21-.6 -^/ror= .06303 06154 
b^ = (.15+4 vC^)/25 = .04519 18359 
b2 = (5.1+11 Vr^)/25 = .31177 75487 
= 1/9 = .11111 11111 

02 - (7+2oV^)/36 = .33052 72081 

03 ={7-Z0^f06)/36 = .05836 16809 

.6- yfM- .35505 10257 OEIGJNAL PA® IS 

b, . -(5.4+19 V^)/25 • -.40216 12204 QUALIIY 

b2 = (20.4+44 V^)/25 = 1.24711 01947 
= 1/9 = .11111 iiri 
32 = (8+5V^)/18 = .51248 58262 
$3 = (8-5 
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9. FOUR FUNCTION INTEGRATION OF x = f{t,x,x) 

This integrator provides fourth order integration for both x and x. 
If X = f(t,x) the integration will be fifth order for both x and x. If 
\ = f^(t) the integrator will be sixth order for both x and k. The inte- 
gration constants are due to Lear. 

/ 

i, • 4Tf(t„.x„,x„) 

l<2 = ATf(tj^ + 62AT, 2^ + ^ K * ^1-1^ 

J <3 = ATf[tj^ + 63 AT, + « 3 AT^ + AT(b^j<^ + 

= ATf[t^ + 64AT, + 6 ^aT^ + AT{c^k^ + C2J^ + 03113 ), ^ 

+ c^Jij + 02^2 

>^+1=4 + AT4 + AT(a^k^ + a2Jl2 + a3il3 + a^l^) 

^+1 ^ * hh ^ 3^3 ■*' ^ 4 ^ 

The integration coefficients are given by 

(5-V^/lO = .27639 32023 
(5+^5)/! 0 = .72360 67977 

1. 

(3-V5)/20 = .038196 60115 
0 

(3+V5)/20 = .26180 33989 
(-l+V5)/4 = .30901 69943 
0 

(3-V5)/4 = .19098 30058 
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= 1/12 = .08333 33333 
ag = (5+V5)/24 = .30150 28324 
03 = (5-V5)/24 = .11516 38343 

«4 = 0 

a.| = (5-V5)/10 = .27639 32023 
= -(5+3V5)/20 = -.58541 01966 
t )2 = (3+Vs)/4 = 1.30901 69944 
^ = -(l-5V?)/4 = 2.54508 49719 
C 2 = -(5+3V?)/4 = -2.92705 09831 
C 3 = (5-Vs)/2 = 1.38196 60113 
= 1/12 = .08333 33333 

^2 = 5/12 = .41666 66667 

63 = 5/12 = .41666 66667 

= 1/12 = .08333 33333 
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10. INTEGRATION OF x = f(t,x,x) + e(t,x,x) 

In many problems there may be a very small perturbing acceleration, 
£(t,x.,A)* adding to a dominate acceleration, f(t,£,^). For example, e^(t) 
could be a small venting acceleration on an orbiting vehicle. £()() could 
be due to higher order harmonics in the earth's gravity field, ^{x) could 
be due to small drag or solar pressure forces acting on a space vehicle. 

When ^ is time consuming to evaluate, it may be desirable to not 
evaluate £ as often as f is evaluated in the integration step. Or, we may 
desire to use an integrator for x = f.(t, 2 ^) when £ has velocity terms in it. 

In these cases we may try the following scheme. 

Integrate f(t,><,£) or £ = f(t,j() normally. Then to the solution add 

4+1 " 4+1 -^V.5’ 4+.5» 4+.5^ 

4+1 " 4+1 -^V.5’ 4+.5» 4+.5^ 

where, we see, £ is evaluated at the midpoint of the integration step. The 
equations for the midpoint are 

t^ + aT/2 

><„ + I - x„) - r - 4> 

2ir ■■ ■ J ■ 4’ ■ 7 ^ 

= (in+1 '*■ 4^/^ AT r 0 or AT small 
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Note that the equation for g is afflicted with roundoff error for small 
aT. a more stable form is 

>^+.5 = Xf, + ^-1 “2-!^ ^ 

in+.5 ~ 4i “2-!^ + •*•) -.25(e^Jc^ + $2^2 ■*■"•) 

where the a^, and are those used in the Nystrom integrator. 
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